#!/usr/bin/perl -w
use strict;
use FindBin;
use lib ("$FindBin::Bin/..", "/net/cpp-group/Leo/bin");
use db_parameters;

#
# connect to genome database
#
my $dbh = connect_to_panda();


print STDERR "\tRead in ortholog_exons\n";

$dbh->do("TRUNCATE orthologs.ortholog_exons;");
$dbh->do("COPY orthologs.ortholog_exons(gene_id1, gene_id2) FROM STDIN");

print STDERR "\t[";
my $cnt_genes = 0;
my $cnt_dots = 0;
while (<>)
{
	
	#gene	ENSCAFG00000000103	ENSG00000074657	2	2	0	2	0	0
	next unless /^gene/;
	++$cnt_genes;
	if ($cnt_genes % 200 == 0)
	{
		++$cnt_dots;
		if ($cnt_dots % 50 == 0)
		{
			print STDERR "\n\t";
			$cnt_dots = 0;
		}
		print STDERR ".";
	}
	my ($gene_id1,  $gene_id2, $selected);
	(undef, $gene_id1,  $gene_id2, undef, undef, $selected) = split /\t/;
	next unless $selected;
	($gene_id1, $gene_id2) = ($gene_id1, $gene_id2) if $gene_id1 gt $gene_id2;
	$dbh->pg_putline("$gene_id1\t$gene_id2\n");
}
$dbh->pg_endcopy();
$dbh->do("VACUUM ANALYSE orthologs.ortholog_exons;");
print STDERR "]\n";


print STDERR "Get list of genes with exons matching its orthologs...\n";
$dbh->do("TRUNCATE orthologs.genes_with_ortholog_exons;");

	my $sql_cmd =<<"PL/SQLCMD";
    INSERT INTO
            orthologs.genes_with_ortholog_exons(gene_id, protocol_id)
        SELECT
                DISTINCT ON (gene_id, protocol_id)
                gene_id, 
                protocol_id
            FROM
            (
                SELECT
                        gene_id1 AS gene_id,
                        protocol_id
                    FROM
                        ortholog_sets os1 JOIN
                        ortholog_sets os2 USING (protocol_id, ortholog_id),
                        ortholog_exons oe
                    WHERE
                        os1.gene_id = oe.gene_id1 AND 
                        os2.gene_id = oe.gene_id2 AND 
                        os1.species <> os2.species
                UNION
                SELECT
                        gene_id2 AS gene_id,
                        protocol_id
                    FROM
                        ortholog_sets os1 JOIN
                        ortholog_sets os2 USING (protocol_id, ortholog_id),
                        ortholog_exons oe
                    WHERE
                        os1.gene_id = oe.gene_id1 AND 
                        os2.gene_id = oe.gene_id2 AND 
                        os1.species <> os2.species
            ) AS a;
PL/SQLCMD
	$dbh->do($sql_cmd);

print STDERR "\tCompleted\n";



$dbh->disconnect();





